Traffic jams and intermittent flows in microfluidic networks 
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We investigate both experimentally and theoretically the traffic of particles flowing in microfluidic 
obstacle networks. We show that the traffic dynamics is a non-linear process: the particle current 
does not scale with the particle density even in the dilute limit where no particle collision occurs. 
We demonstrate that this non-linear behavior stems from long range hydrodynamic interactions. 
Importantly, we also establish that there exists a maximal current above which no stationary particle 
flow can be sustained. For higher current values, intermittent traffic jams form thereby inducing the 
ejection of the particles from the initial path and the subsequent invasion of the network. Eventually, 
we put our findings in the broader context of the transport proccesses of driven particles in low 
dimension. 
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Hundreds of industrial and biological processes ulti- 
mately rely on the transport of particle or droplet sus- 
pensions in channels or obstacle networks. Prominent ex- 
amples include: particle filtration [3 [2] , size separation 
of colloids and polymers pi l^, emulsion flows through 
porous rocks during enhanced oil recovery, droplet trans- 
port in microfluidic chips [5J [5] and, obviously, blood mi- 
croflows [7|. For all these examples, the particle size typ- 
ically compares with the distance between the obstacles, 
or walls, forming the network. Therefore, even in the 
dilute limit, the particles locally reduce the network con- 
ductivity, thereby inducing a dynamic redistribution of 
the fluid flow. This non-trivial interplay between the po- 
sition of the advected particles and the driving flow fleld 
makes the description of the particle traffic a challenging 
task. Whereas all the above examples typically involve 
very large networks, most of the recent investigations of 
such "dynamic" networks have been hitherto focused on 
droplet traffic in minimal microfluidic set-ups composed 
of a single loop or of a single T-junction. Nonetheless, 
despite the apparent simplicity of these geometries, a 
rich variety of periodic and aperiodic drop dynamics have 
been reported due to finite time anti-correlations of the 
droplets trajectories at the junctions [8HTi]. 

In this letter, we aim at providing a generic descrip- 
tion of the trafficking dynamics in large fluidic networks. 
More precisely, we focus on a wide class of prototypal 
microfluidic networks where a single particle would fol- 
low a unique deterministic trajectory connecting a single 
inlet to a single outlet, see Fig. [l^. Increasing the par- 
ticle injection rate, j, we observed that above a critical 
value j* , a traffic jam builds up at the inlet. This lo- 
cal increase of the density results in the ejection of the 
particles, which then invade a larger fraction of the net- 
work and explore alternative exit paths, Fig.jTJj, c, and d. 
Surprisingly, this invasion transition occurs at a very low 
apparent volume fraction for which the particles do not 
experience any collision. In the high current regime, we 
also report on the collective ejection of the particles out 



of the initial lane which results in an intermittent traf- 
fic flow. We combine experimental investigations with a 
minimal theoretical model to uncover the physical origin 
of the jamming and of the invasion dynamics of fluidic 
networks. 

Our experiment consists in flowing periodic sequences 
of droplets through an obstacle network inside a microflu- 
idic channel. The network is made of a 15x33 square lat- 
tice of cylindrical posts. Figure [l] The diameter of the 
posts is 150 /im, the lattice spacing, (., is 150 /im and the 
channel height is 80 //m. Changing the lattice size and 
the post shape did not yield any qualitatively different 
results (results not shown). We chose this design so that 
a single droplet follows a unique deterministic trajectory 
along the central lane, which is the shortest path con- 
necting the entrance to the exit. We used water droplets 
(viscosity r/ =1 mPa.s) as the advected particles rather 
than solid microbeads to avoid irreversible clogging of 
the network. The droplets were made and transported 
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FIG. 1: Water droplets in hexadecane oil flowing in a square 
lattice of obstacles, (a) Imposed droplet current: j = 0.4 Hz. 
(b) Droplet current: j = 0.7 Hz. A jam forms at the entrance 
of the network, (c) Droplet current: j = 2 Hz. (d) Very 
high current phase. The droplets invade the whole network. 
Length of the network: 9.75 mm. 
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FIG. 2: Number of droplets in the central lane: A'^ (•), and 
number of droplets in the other lateral lanes, A'', (D) plotted 
versus j (experiments). Full red line: theoretical prediction 
for the variations of A^ according to Eq. Is] (best fit). The 
dashed vertical line indicates where the invasion transition 
occurs: j = j* ■ 



FIG. 3: Time averaged density profiles normalized by the 
close packing density l/(2a) in the central lane for different 
imposed droplet currents j. The position is normalized by the 
lattice spacing, £. 



in a mixture of hexadecane (viscosity i] — 2inPa.s) and 
span 80 surfactant (3wt%). To produce the droplets in 
situ we used a home made drop-on-demand device in- 
troduced previously in [15; . This device allowed us to 
control independently and accurately: the droplet size, 
the flow rate of the continuous phase and the droplet 
current: j, i.e. the number of droplets per second in- 
jected in the device, j was varied over 2 decades from 
0.05 to 5 Hz. In order to ensure a long term stability 
and a fast dynamic response of the droplet emitter, we 
used a microfluidic sticker made of a photocurable op- 
tical adhesive (N0A81, Norland Products) for the main 
microfluidic channel |16j . In all the experiments, we used 
droplets with a radius a = 89 ± 4 /xm and the oil flow rate 
was set so that the velocity of a single droplet in the cen- 
ter of the channel was v'^ = 350 /im/s. Thus, both the 
Reynolds number, Re ~ 10~^ and the capillary number, 
Ca ^ 10~^ were fixed at a constant small value in all our 
experiments. 

We first show in Figure 2, how the numbers of droplets 
in the central lane, N, and outside the central lane, N, 
change as the imposed droplet current j increases. We 
also plot in Fig. [3] the density field in the central lane, 
p{x), normalized by its maximal close-packing value, 
l/(2a), where again a is the droplet diameter. p{x), N 
and N, are time-averaged quantities systematically mea- 
sured over 25-minute experiments, performed 10 minutes 
after the injections of the first droplets at a given current 
j- 

In the very small current limit, j < 0.3 Hz, N{j) in- 
creases linearly with j, Fig. 2. The particles keep on fol- 
lowing the central lane as an isolated particle would do. 
The local density p{x) simply reflects the x component of 
the velocity field, Vp{x), of the continuous phase in a par- 



ticle free network, Fig. [3J The imposed droplet current, 
and the local density are linearly related: j oc p{x)vp{x) . 
The variations of p close to the entrance and the exit are 
due to the compression (resp. elongation) of the bare 
fiow field at the entrance (resp. exit) of the network. In 
all that follows, we ignore this linear hydrodynamic ef- 
fect and focus only on the bulk properties of the p field. 
Increasing the current, the traffic flow crosses-over to a 
qualitatively different regime. 

For 0.3 Hz< j < 0.6 Hz, the droplets stiU follow 
the central lane and p{x) retains the same fiat profile. 
Fig. [SJ However, the N{j) curve significantly deviates 
from its initial linear behavior and its variations diverge 
as the current approaches a well-defined "critical value" 
j* — 0.6 Hz, see Fig. [2J In fact, these non-linear varia- 
tions reflect the slowing down of the droplets as their lin- 
ear density increases. We now stress on one of our most 
important result: j* is the maximal stationary particle 
current that can be sustained in the network. Indeed, as 
j reaches j* , the uniform and stationary droplet distribu- 
tion becomes unstable: a jam forms a the entrance of the 
network, the droplets escape the central lane and invade 
the network. Figs. [T] and [2j Surprisingly, this instabil- 
ity takes place though the droplet density remains much 
smaller than the close-packing density, p{x) < 0.7/(2a), 
and no droplet collision occurs. We can thus unambigu- 
ously infer that long range hydrodynamic interactions 
are shaping the non-linear (j, N) constitutive relation in 
Fig. [2] and are responsible for the unexpected destabiliza- 
tion of the traffic flow at j*. Before describing the higher 
current regimes, we now introduce a theoretical model to 
uncover the origin of these effective hydrodynamic inter- 
actions. 

Since we investigate rather dilute systems, we use a 
continuous description that substantially simplifies the 
algebra. This simplification is done without loss of gen- 



erality in the far field limit. Firstly, we model the ob- 
stacle network by a lattice of ID channels solely charac- 
terized by their local hydrodynamic conductivity G. All 
the specifics of the obstacle shape are hidden in this pa- 
rameter which relates the local pressure gradient to the 
local fluid velocity: VF(r) = — G(r)VP(r), where r is the 
position of the considered elementary channel in the 2D 
network. Using the incompressibility relation, V-vp = 0, 
we then obtain the equation for the pressure field: 



V • [G(r)VP(r)] = 0, 



(1) 



with the boundary condition: VP = —{v-p/G)x, where 
vf is the constant fluid velocity far from the particles. For 
a regular lattice of obstacles, G is constant. Secondly, 
following Jousse and coworkers [51 [TO] we describe the 
advected droplets as constant pointwise perturbations to 
the local conductance. Indeed, the droplet size is smaller 
than the lattice spacing, £, which defines the cut-off of 
our continuous theory. More precisely, for N particles 
located at rj(i) = {xi,0): 



G{r) = G 



N 



l_52^5(r-ri(t)) 



(2) 



where the S function is regularized at short distance: 
b^S{0) ^ {b/f-Y . Here b/i measures how much a droplet 
locally hinders the fiuid flow. Thirdly, the last ingredi- 
ent of this model is the droplet advection rule. For sake 
of simplicity, we assume that a droplet located at Yi{t) 
moves with a velocity r^ = ^wp{vi) where /x is a constant 
mobility coefficient. As the flow is potential, we can now 
take advantage of an obvious electrostatic analogy. In- 
deed, for iV = I, Eqs.[T] and [2] correspond to the equation 
for the electric potential induced by a pointwise dielectric 
particle polarized by a homogeneous electric field [T7] . 
Therefore, we readily deduce that, in the far field limit, 
a single droplet acts as a source dipole oriented along 
the X-axis. More quantitatively, we can compute the am- 
plitude of the corresponding fiow disturbance along x at 
r = (a;, 0) by solving perturbatively Eq. [Tj At first or- 
der in 6^: Wdip(x;ii) = ii6^(I - ^)/{2T:^)dxj^^\. We 
stress on the sign of this perturbation, which is opposed 
to the driving flow along the x-axis. Therefore, it re- 
duces the velocity of two approaching droplets whatever 
their relative position. We now write down an effective 
equation for the dynamics of N aligned particles. To 
do so we exploit the pairwise additivity of the hydro- 
dynamic interactions within the above dipolar approxi- 
mation. Hence, the equations of motion take the simple 
form: dtXi{t) = Vy) + ^jLY,j^i^Aip{xj-Xi;Xj), where Vj^ is 
the velocity of a single advected droplet. These N equa- 
tions can be solved analytically when the droplets are 
equally separated along a single lane by a distance L/N, 
where L is the system size, and all moving at the same 
speed, Xi = vd. This corresponds to the only stationary 



solution we observed in our experiments. Summing over 
all the dipolar interactions is straightforward and yields 
a non-linear relation between the density p — N/L and 
the uniform particle current: j = pwD- In the large N 
limit: 



J 



P^D 



I+f(I-A*)&2p2 



(3) 



As we measured the velocity w^, independently , we are 
left with a single fitting parameter. For (1 — /i)6^ ~ 
0.17 mm~^, our first order prediction is in excellent agree- 
ment with our experimental data, see the red curve in 
Fig. [2| Our second central theoretical result is that the 
the j(p) function defined in Eq. 3 has an absolute max- 
imum at J = j* — 0.65 Hz, implying that there cannot 
exist any stationary solution for j > j* ■ This result pro- 
vides a quantitative explanation for the traffic destabi- 
lization observed experimentally above this specific cur- 
rent value. 

We shall add that the above model is not restricted to 
2D square lattices. It is indeed straightforward to extend 
it to any homogeneous potential Stokes flow, which en- 
compasses flows in bidimensional and tridimensional reg- 
ular networks of arbitrary symmetry and droplet trans- 
port in Hele-Shaw geometries. For instance, our predic- 
tion is in agreement and generalize the hydrodynamic 
slowing down of droplets in obstacle-free shallow chan- 
nels reported in [TOl [T^ . 

In the last part of this letter, we focus on the traffic 
dynamics in the high current regime j > j* ■ For current 
values larger but close to j*, typically for j* < j < I Rz, 
a jam forms at the entrance of the microfiuidic network: 
the instantaneous density field p{x, t) becomes asymmet- 
ric and its maximum value continuously increases until 
the droplets locally contact each other. Figs, [l] and [3J 
Subsequently, droplets are ejected toward the upper and 
lower lanes and start invading the network. Figs, [l] and [2] 
Consequently, the local density in the central lane is re- 
duced but it does not reach a stationary value, as the 
jam reforms slowly and disrupts again quickly: the traf- 
fic flow becomes intermittent. To gain a better insight on 
this invasion dynamics, we plot in Fig. [4] a, the variations 
of the instantaneous number of droplets in the central 
lane, N{t). We clearly distinguish two time scales in this 
asymmetric signal: the slow one corresponds to the jam 
formation, where N(t) increases, and the faster one cor- 
responds to the time during which droplets are coopera- 
tively ejected from their initial path and N{t) decreases. 
In addition, we plot in Fig. [4|b, the probability distribu- 
tion, P{N), of the number of particles traveling outside 
the central lane for j = 0.85 Hz. This distribution is 
clearly non-Gaussian, which is another signature of the 
intermittent dynamics and of the correlated ejection of 
groups of droplets. Rare events, during which groups of 
particles travel cooperatively through the lateral lanes, 
broaden the distribution and increase the mean number 




FIG. 4: (a) Variations of N{t) in the central lane for j = 0.85 
Hz. The red arrow indicates cooperative-ejection events. (b) 
Histogram of the number of droplets flowing outside the cen- 
tral lane P{N), computed for 1050 droplets injected in the 
network at a rate j — 0.85 Hz. Vertical line: mean value (c) 
N{t) curve for j = 1.2 Hz, same time interval as in (a), (d) 
Histogram P{N), computed for 3000 droplets injected in the 
network at a rate j = 1.2 Hz. Sketch of the effective hydro- 
dynamic interactions between 3 advected droplets flowing on 
two adjacent lanes. 



of particles ejected out of the central lane above its typi- 
cal value. The cooperative ejection of the droplets can be 
qualitatively understood by looking at the symmetry of 
the dipolar flow field. As sketched in Fig. |4j the droplets 
moving above and below the central lane increase the 
transverse flow downstream, thereby inducing the pump- 
ing of the droplets from the central to the side lanes, see 
the picture in FigHt. In turn, this pumping mechanism 
is amplified as the number of droplets on side lanes in- 
creases, which results in an unstable dynamics yielding 
the correlated ejection of several droplets when the local 
density reaches a sufhciently high value. 

Finally we briefly describe the high current regime: 
j > 1 Hz. For imposed currents much larger than j* , the 
traffic flow remains non-stationary. However, the jams 
forms and disrupts continuously without any clear sepa- 
ration between the time scale of these two processes, see 
Fig. |4| c. Consequently, the distribution of N becomes 
more symmetric, and seems to converge toward a Gaus- 
sian behavior in the large j limit, Fig. |4]d. In addition, 
we found that above j* , the system self-organizes so that 
the density field along the initially preferred path con- 
verges and becomes independent of the injected current 



value. To go beyond this first attempt at understanding 
the traffic dynamics above j*, a more thorough statisti- 
cal analysis would be required. Work along this line is in 
progress. 

In summary, we have demonstrated that when particles 
are advected in a fluidic network, they locally reduce the 
network conductance thereby inducing non-local dipolar 
interactions between the moving objects. The main con- 
sequence of this hydrodynamic coupling is the limitation 
of the maximal stationary current, j* that can be sus- 
tained in the network. Above j*, the initially preferred 
path becomes unstable and the particles intermittently 
invade the whole network. To close this letter, we note 
that the destabilization and the jamming of stationary 
traffic flows go beyond the present study. The exis- 
tence of a maximal current is indeed a very generic fea- 
ture reported for numerous classes of theoretical models 
of driven interacting-particle systems, including discrete 
exclusion processes |20] and phenomenological continu- 
ous models of self-driven particles (pedestrian or vehi- 
cles) [3T]. All these models possess a common feature: a 
non-monotonic constitutive relation relating the current 
to the particle density . To our knowledge, we provide the 
first experimental evidence of such a maximal current. 

We thank Julien Tailleur for insightful comments and 
suggestions. Raphael Voituriez is also ackowledged for 
an illuminating discussion. 
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